Vatmal block length
data
We’ve taken all non-road landuse polygons in vatmal (20622 polygons), and dissolved the נorders in blocks(5606 blocks). this is an example from arara:
analysis
we have few ways to take blocks into accounts: 1) as they are 2) as Minimum Oriented Rectangle(MOR), which is very similar to a bounding box, just with an angle 3) as Convex Hull(CH) - a shape that “wrap” the block as a present, and makes sure it does not have any reflex angles
it is imperative to identify the right scope of analysis. below are graphs that show the cumulative distribution of relative growth in block size in CH and MOR:
blocks which were less than 2000 square meters were removed from the analysis.
calculate_segment_angles_directions <- function(geometry) {
if (st_is(geometry, "LINESTRING")) {
# Extract the vertices from the linestring
vertices <- st_coordinates(geometry)
# Calculate the differences between consecutive vertices
diffs <- diff(vertices)
# Calculate the angles between consecutive segments using atan2
angles <- atan2(diffs[c(nrow(diffs),1:(nrow(diffs)-1)), "Y"], diffs[c(nrow(diffs),1:(nrow(diffs)-1)), "X"]) - atan2(diffs[, "Y"], diffs[, "X"])
angles <- angles * 180 / pi
# Adjust the angles to be between 0 and 360 degrees
angles <- ifelse(angles < 0, angles + 360, angles)
# Calculate the lengths of segments before and after each angle
lengths_after <- sqrt(diffs[, "X"]^2 + diffs[, "Y"]^2)
lengths_before <- sqrt(diffs[c(nrow(diffs),1:(nrow(diffs)-1)), "X"]^2 + diffs[c(nrow(diffs),1:(nrow(diffs)-1)), "Y"]^2)
# Calculate the directions based on the angles
directions <- ifelse(angles > 180, "Left","Right")
angles <- ifelse(angles > 180, 360 - angles, angles)
# Return a data frame with the angles, lengths before, lengths after, and directions
return(data.frame(angles = angles,
lengths_before = lengths_before,
lengths_after = lengths_after,
x= vertices[-nrow(vertices),1],
y= vertices[-nrow(vertices),2],
directions = directions))
}
if (st_is(geometry, "POLYGON")) {
# Extract the exterior ring of the polygon
ring <- st_cast(geometry, "LINESTRING")
# Calculate the angles, lengths before, lengths after, and directions using the linestring calculation
result <- calculate_segment_angles_directions(ring)
# Return the result
return(result)
}
# Return NULL if the geometry type is not supported
return(NULL)
}we are splitting the blocks perimeter by left and right turns into components. there were parameters involved: if the length of a component was less than 26 meters and the angle turn sum was less than 40 degrees, it was united with the component beforehand afterwards, we carried out dbscan to cluster nearby points (less than 5 meters) in the block, to better understand where cuts are necessary eventually, we count the number of components in each block. we have 4326 blocks with up to 2 components. and 284 with 3 or more. their analysis is harder. however we will try to analyze all according to all of the methods.
.
1 2 4 6 8 10 14 16 40 42
3210 1116 216 46 15 2 2 1 1 1 
# download_this(ggpols2_res %>% bind_cols(pols2) %>% mutate(geometry = map_chr(geometry,st_as_text,digits = 8)),
# output_name = "pols2",
# output_extension = ".csv",
# button_label = "Download polygons analysis",
# button_type = "warning",
# has_icon = TRUE,
# icon = "fa fa-save",
# id = 'pols' # add an id to this button
# )
ggpols2_res %>% bind_cols(pols2) %>% st_sf() %>% st_write("pols.gpkg",delete_layer =T)Deleting layer `pols' using driver `GPKG'
Writing layer `pols' to data source `pols.gpkg' using driver `GPKG'
Writing 4610 features with 53 fields and geometry type Polygon.
# download_this(ggmor2_res %>% bind_cols(pols2) %>% mutate(geometry = map_chr(geometry,st_as_text,digits = 8)),
# output_name = "mor2",
# output_extension = ".csv",
# button_label = "Download MOR analysis",
# button_type = "warning",
# has_icon = TRUE,
# icon = "fa fa-save",
# id = 'mor' # add an id to this button
# )
ggmor2_res %>% bind_cols(mor2) %>% st_sf() %>% st_write("mor.gpkg",delete_layer =T)Deleting layer `mor' using driver `GPKG'
Writing layer `mor' to data source `mor.gpkg' using driver `GPKG'
Writing 4610 features with 11 fields and geometry type Polygon.
# download_this(ggch2_res %>% bind_cols(pols2) %>% mutate(geometry = map_chr(geometry,st_as_text,digits = 8)),
# output_name = "ch2",
# output_extension = ".csv",
# button_label = "Download CH analysis",
# button_type = "warning",
# has_icon = TRUE,
# icon = "fa fa-save",
# id = 'ch' # add an id to this button
# )
ggch2_res %>% bind_cols(ch2) %>% st_sf() %>% st_write("ch.gpkg",delete_layer =T)Deleting layer `ch' using driver `GPKG'
Writing layer `ch' to data source `ch.gpkg' using driver `GPKG'
Writing 4610 features with 12 fields and geometry type Polygon.